Introduction: The Number We’ve Seen Since Part 2
Every time you look at a PHREEQC output, you’ve probably wondered about the lengthy list of minerals and numbers at the end. Up until now, we might have turned a blind eye to this specific section. Let’s finally decode it:
Calcite 0.39 -7.83 -8.22
Dolomite 0.65 -16.65 -17.09
Gypsum -0.48 -5.28 -4.80
Kaolinite 2.64 -8.11 -8.11
Quartz 0.18 -3.88 -4.00
Pyrite -2.34 -30.44 -28.10
This is the SI (Saturation Index)—arguably the most frequently referenced number in a PHREEQC output. If you can read it correctly, you can tell at a glance “which minerals the water is dissolving” and “which minerals the water is trying to precipitate.”
- The definition of SI and the meaning of IAP and \(K_{sp}\)
- An intuitive understanding of what \(SI > 0\), \(SI = 0\), and \(SI < 0\) signify
- Deciphering SI patterns in natural waters (acid rain, groundwater, hot springs, seawater)
- How to extract the SI for multiple minerals simultaneously using PHREEQC code
- Why Kaolinite SI often shows high values and its geochemical implications
- Dealing with uncertainties in SI (measurement errors and temperature corrections)
- Practical applications for water treatment and scale prediction
Theory: Definition of SI
Ion Activity Product (IAP) and Solubility Product (\(K_{sp}\))
Let’s take the dissolution equilibrium of a mineral as an example. The dissolution of Calcite (CaCO₃) is written as:
\[\text{CaCO}_3 \text{(s)} \rightleftharpoons \text{Ca}^{2+} + \text{CO}_3^{2-}\]
The equilibrium constant (Solubility Product) is:
\[K_{sp} = a_{\text{Ca}^{2+}} \cdot a_{\text{CO}_3^{2-}} = 10^{-8.48} \quad (25°C)\]
The Ion Activity Product (IAP) in the actual solution is:
\[\text{IAP} = a_{\text{Ca}^{2+}}^{\text{actual}} \cdot a_{\text{CO}_3^{2-}}^{\text{actual}}\]
Definition of SI
\[SI = \log_{10}\left(\frac{\text{IAP}}{K_{sp}}\right) = \log(\text{IAP}) - \log(K_{sp})\]
Mineral will dissolve
(if present)
Acid rain, Soft water
No net dissolution
or precipitation
Limestone aquifer
Mineral will precipitate
(if nucleation occurs)
Hot springs, Brines
An \(SI > 0\) does not guarantee that a mineral will precipitate instantly. Nucleation requires activation energy, and the precipitation of inorganic minerals can sometimes take hours, years, or even millennia (a state of sustained supersaturation). Conversely, even with an \(SI < 0\), dissolution is not instantaneous and depends on kinetics. The SI merely indicates the thermodynamic driving force.
Visualizing the Concept of SI
Figure 1. The SI Concept — An analogy between Water Level (IAP) and a Dam (Ksp)
Calculating Calcite SI by Hand
To verify the theory, let’s manually calculate the SI of Calcite in seawater.
# From Part 9, γ(Ca²⁺) = 0.28 in seawater
m(Ca²⁺) = 0.0103 mol/kg
a(Ca²⁺) = 0.28 × 0.0103 = 0.00288
# Activity of CO₃²⁻ (From PHREEQC output)
a(CO₃²⁻) = 2.37×10⁻⁵ (pH 8.22, Alk 2.3 mM)
IAP = a(Ca²⁺) × a(CO₃²⁻)
= 0.00288 × 2.37×10⁻⁵
= 6.83×10⁻⁸
log IAP = log(6.83×10⁻⁸) = -7.17
log Ksp (Calcite, 25°C) = -8.48
SI = -7.17 - (-8.48) = +0.72 ← Matches PHREEQC exact output of +0.721
This means that seawater is thermodynamically primed to precipitate calcite. However, it does not spontaneously turn into a milky solid mass because of the activation energy for nucleation and the inhibiting effect of Mg²⁺ in seawater.
PHREEQC Code (Complete Workflow)
Setup: Batch Calculating SIs for 5 Natural Water Types
To correctly calculate the SI for Calcite or Kaolinite, you must define Al and Fe in all solutions (even if only in trace amounts). If they are strictly undefined (concentration = 0), PHREEQC cannot compute the IAP and will return an SI = -999.999 (incalculable). Additionally, the Carbonate system (Calcite) requires either C(4) or Alkalinity to be defined.
# ============================================================
# DeepFlow Part 10 - Batch SI Calculation
# 5 Natural Waters × Calcite, Quartz, Kaolinite, Goethite
# ============================================================
SOLUTION 1 "Acid Rain"
temp 15
pH 4.5
pe 12
units mol/kgw
S(6) 1e-6
N(5) 5e-6
Ca 2e-6 # Trace Ca from atmosphere
C(4) 1e-5 # Atmospheric CO2 (pCO2 ≈ 400ppm equivalent)
Si 1e-5 # Trace Si from atmosphere
Al 1e-7 # Trace Al before weathering
Fe(2) 1e-7
Fe(3) 1e-10
SOLUTION 2 "Limestone Groundwater"
temp 12
pH 7.4
pe 4
units mol/kgw
Ca 2.1e-3
Mg 0.8e-3
Na 0.5e-3
K 0.1e-3
Cl 0.4e-3
Alkalinity 5.2e-3 as HCO3
S(6) 0.3e-3
Si 0.1e-3
Al 1e-7
Fe(2) 1e-7
Fe(3) 1e-10
SOLUTION 3 "Granite Groundwater"
temp 10
pH 6.2
pe 6
units mol/kgw
Ca 0.4e-3
Mg 0.1e-3
Na 1.2e-3
K 0.3e-3
Cl 0.8e-3
Alkalinity 1.5e-3 as HCO3
Si 0.3e-3
Al 2e-6 # Sourced from feldspar weathering
S(6) 1e-5
Fe(2) 1e-7
Fe(3) 1e-10
SOLUTION 4 "Carbonate Hot Spring"
temp 50
pH 6.8
pe 2
units mol/kgw
Ca 5.0e-3
Mg 1.0e-3
Na 8.0e-3
K 0.5e-3
Cl 2.0e-3
Alkalinity 12.0e-3 as HCO3
S(6) 1.0e-3
Fe(2) 0.5e-3
Si 1.0e-3
Al 5e-7
SOLUTION 5 "Seawater"
temp 25
pH 8.22
pe 8.451
units mol/kgw
Na 0.4689
Mg 0.0528
Ca 0.0103
K 0.0102
Cl 0.5453
Alkalinity 2.3e-3 as HCO3
S(6) 0.0283
Si 1.0e-4
Al 1e-8
Fe(2) 1e-8
Fe(3) 1e-10
# -------- SI Calculation & Output --------
SELECTED_OUTPUT 1
-file si_comparison.txt
-solution true
-pH true
-temperature true
-ionic_strength true
-saturation_indices Calcite Quartz Chalcedony Kaolinite Goethite
USER_PUNCH 1
-headings Solution Temp pH Ionic_Strength SI_Calcite SI_Quartz SI_Chalcedony SI_Kaolinite SI_Goethite
-start
10 PUNCH DESCRIPTION, TEMP, -LA("H+"), MU, \
SI("Calcite"), SI("Quartz"), SI("Chalcedony"), SI("Kaolinite"), SI("Goethite")
-end
USE solution 1
END
USE solution 2
END
USE solution 3
END
USE solution 4
END
USE solution 5
END
Results and Interpretation
The table below shows the execution results from the script above.
| Solution | pH | SI Calcite |
SI Quartz |
SI Kaolinite |
Diagnosis |
|---|---|---|---|---|---|
| Acid Rain | 4.5 | −10.13 | −0.87 | −5.52 | Strongly dissolves all minerals. Carbonates and Kaolinite are undersaturated. |
| Limestone GW | 7.4 | +0.15 | +0.18 | +2.64 | Calcite & Quartz are near equilibrium. Kaolinite is supersaturated. |
| Granite GW | 6.2 | −2.23 | +0.69 | +6.96 | Carbonate undersaturated. Quartz supersat. Extreme supersat for Kaolinite (see below). |
| Hot Spring | 6.8 | +0.65 | +0.64 | +2.76 | Calcite & Quartz supersaturated. High risk for scale formation (SI > 0.5). |
| Seawater | 8.22 | +0.72 | +0.04 | −2.77 | Calcite supersat, Quartz near equilibrium. Kaolinite undersaturated (due to high ionic strength). |
Why is the Kaolinite SI So High?
The extremely high value of +6.96 for Kaolinite in Granite Groundwater might seem like an error, but it has profound geochemical significance.
The dissolution reaction of Kaolinite is:
\[\text{Al}_2\text{Si}_2\text{O}_5(\text{OH})_4 + 6\text{H}^+ \rightleftharpoons 2\text{Al}^{3+} + 2\text{H}_4\text{SiO}_4 + \text{H}_2\text{O}\]
The \(\log K\) for this reaction is approximately \(-7.4\) (in phreeqc.dat). The IAP is determined by the activities of Al³⁺, Si, and H⁺. Because the contribution of H⁺ (pH) is dominant, even a trace amount of Al (e.g., 2 μmol/L) at pH 6.2 will result in a massive IAP.
An SI of around +7 is completely natural in the weathering front of granitic terrain. It signifies that the alteration of feldspar into kaolinite is strongly driven thermodynamically. In reality, kaolinite forms as a secondary mineral on the surface of dissolving feldspars. Although it is thermodynamically favored, its formation is constrained by the slow kinetics of feldspar dissolution and nucleation, making it a mineral that forms over geological timescales.
Application 1: Scale Prediction (Calcite Precipitation)
In geothermal power plants and hot springs, as water cools, the SI changes, leading to potential scale formation. Our modeled Hot Spring water is already supersaturated at its 50°C source (SI = +0.65). We can track how this changes as it cools using phreeqpy (Python binding).
import os
import pandas as pd
from phreeqpy.iphreeqc.phreeqc_dll import IPhreeqc
import matplotlib.pyplot as plt
def run_scale_prediction():
current_dir = os.path.dirname(os.path.abspath(__file__))
db_path = os.path.join(current_dir, "phreeqc.dat")
phreeqc = IPhreeqc()
phreeqc.load_database(db_path)
input_script = """
SOLUTION 1 "Hot Spring Source"
temp 50
pH 6.8
units mol/kgw
Ca 5.0e-3
Mg 1.0e-3
Na 8.0e-3
K 0.5e-3
Cl 2.0e-3
Alkalinity 12.0e-3 as HCO3
S(6) 1.0e-3
Fe(2) 0.5e-3
Si 1.0e-3
Al 5e-7
REACTION_TEMPERATURE 1
50 30 15
SELECTED_OUTPUT 1
-reset false
-temperature true
-ph true
-si Calcite Aragonite
"""
phreeqc.run_string(input_script)
output = phreeqc.get_selected_output_array()
if output:
return pd.DataFrame(output[1:], columns=output[0])
return "No results."
if __name__ == "__main__":
results = run_scale_prediction()
print(results)
plt.figure(figsize=(8, 5))
plt.plot(results['temp(C)'], results['si_Calcite'], marker='o', label='Calcite')
plt.plot(results['temp(C)'], results['si_Aragonite'], marker='s', label='Aragonite')
plt.axhline(0, color='red', linestyle='--', alpha=0.5)
plt.title("Temperature vs Saturation Index (SI)")
plt.xlabel("Temperature (°C)")
plt.ylabel("Saturation Index (SI)")
plt.gca().invert_xaxis()
plt.legend()
plt.grid(True)
plt.show()pip install phreeqpyEnsure phreeqc.dat is in the same directory as your script. Avoid special characters or spaces in the file path.
Figure 2. Decrease in SI(Calcite) as Hot Spring Water Cools — Even as the SI decreases, it remains in the supersaturated zone, meaning scaling risk persists.
Calcite is known for its retrograde solubility, meaning it usually becomes less soluble as temperature increases, and more soluble as it cools. However, in our Python code, we only provided an initial Alkalinity of 12 mmol/L as the carbon source. We did not add excess C(4) (like CO₂ gas from the atmosphere).
Under these closed conditions, as the water cools, CO₂ becomes more soluble (Henry’s Law), causing the balance to shift towards aqueous \(\text{H}_2\text{CO}_3^*\). This causes \(\text{CO}_3^{2-}\) to shift toward the \(\text{HCO}_3^-\) side (association) and decrease in concentration. The resulting drop in IAP leads to the SI decreasing alongside temperature (which looks like normal solubility).
Retrograde behavior becomes strongly apparent in open systems with fixed \(p\text{CO}_2\) (e.g., exposed to the atmosphere) or in closed systems with very high initial total carbon. Accurate scale diagnosis requires precise knowledge of the field’s \(p\text{CO}_2\) conditions!
Application 2: Diagnosing Water Types by “Reading” SI Patterns
By looking at the pattern of Saturation Indices, you can deduce the geological history of the water.
SI(Quartz) = +0.18 (Near Eq)
SI(Kaolinite) = +2.64 (Supersat)
Ca–HCO₃ type water. It is close to equilibrium with calcite. The longer the residence time, the closer SI(Calcite) gets to 0. Kaolinite supersaturation indicates ongoing feldspar weathering into clay.
SI(Quartz) = +0.69 (Supersat)
SI(Kaolinite) = +6.96 (Extreme Supersat)
Na–HCO₃ type. Weathering of feldspar supplies Na, Si, and Al. The massive SI for Kaolinite is a hallmark of active feldspar weathering and clay formation.
SI(Quartz) = +0.64 (Supersat)
SI(Kaolinite) = +2.76 (Strong Supersat)
High alkalinity (12 mmol/L) drives carbonate supersaturation. Extreme caution is needed for calcite scaling in pipes and heat exchangers.
SI(Quartz) = +0.04 (Near Eq)
SI(Kaolinite) = −2.77 (Undersat)
Carbonates are supersaturated but rarely precipitate spontaneously without biological intervention. Due to high ionic strength, Al activity is suppressed, keeping Kaolinite undersaturated.
Uncertainties: Handling Measurement Errors
A field measurement error of ±0.1 in pH translates directly to roughly ±0.1 uncertainty in SI. A temperature error of ±2°C causes an SI variance of ±0.1 to 0.2 in carbonate systems. Therefore, practically speaking, any SI between -0.2 and +0.2 should be interpreted as “essentially at equilibrium.”
| SI Range | Geochemical Meaning | Practical Interpretation |
|---|---|---|
| SI < -1 | Strongly Undersaturated | Fast dissolution for carbonates / Kinetic limits for silicates |
| -1 ≤ SI < -0.2 | Undersaturated | No scaling; slow long-term mineral depletion |
| -0.2 ≤ SI ≤ +0.2 | Near Equilibrium | The mineral is likely controlling the system chemistry |
| +0.2 < SI ≤ +1 | Mildly Supersaturated | Potential precipitation (scaling starts for carbonates) |
| SI > +1 | Strongly Supersaturated | High scaling risk / Extreme thermodynamic drive (though kinetics may still limit actual precipitation) |
References
Other articles in this series:
Using the REACTION block, we will incrementally dissolve rock fragments and track how pH, SI, and species concentrations evolve. It’s a simulation of “rainwater contacting limestone and transforming into groundwater.” This will integrate all the knowledge we’ve built up from Part 3 through 10.
Other articles in this series:
- #1 Installation and Initial Calculation
- #2 Analyzing Seawater with Speciation
- #3 Mineral Equilibrium and Temperature Effects
- #4 Calcite–CO₂ Interaction (Open vs. Closed Systems)
- #5 Mixing Groundwater and Seawater
- #6 Pyrite Oxidation and AMD Formation
- #7 Solubility Diagrams (Gibbsite)
- #8 Visualization with Python
- #9 Ionic Strength and Activity Coefficients
- #10 Mastering the Saturation Index (SI) (This article)
- #11 Reaction Path Modeling
DeepFlow | Science beneath the surface